      Subroutine GradCorrM1(IB, M)  ! 2D version only
      include "param.fi"
      include "common.fi"
      INTEGER IB
      INTEGER IJK, II, JJ, KK, DELI, DELJ, DELK, I, J, K, IP, IB2
      REAL TMPX, TMPY, TMPZ, KX, KY, KZ, Distance
      REAL M(2,2)

      Do I = 1, 2
        Do J = 1, 2
          M(I, J) = 0.0
        Enddo
      Enddo

        ! Get Hash Index
        IJK = HASS(IB)
        Call HASHINDEX13(IJK, II, JJ, KK)
        Do DELI = -1, 1
        I = HASHIND(II + DELI)
        Do DELJ = -1, 1
        If (DIMM .EQ. DIMM2D .AND. DELJ .NE. 0) Cycle
        J = HASHIND(JJ + DELJ)
        Do DELK = -1, 1
        K = HASHIND(KK + DELK)
          Do IP = 1, HASSNUMS(I,J,K)
            IB2 = HASSIND(IP, I, J, K)
            Call Distance2P_VEC(X1(IB2),X2(IB2),X3(IB2),
     &                          X1(IB),X2(IB),X3(IB),
     &                          TMPX, TMPY, TMPZ)
            Distance = SQRT(TMPX*TMPX+ TMPY*TMPY+ TMPZ*TMPZ)
            If (Distance.LE.HSUP .AND. Distance .GT. 1E-10 
     &          .AND. IB.NE.IB2) Then
                Call KernelGradient(TMPX, TMPY, TMPZ, KX, KY, KZ)
                M(1,1) = M(1,1) - KX * Weight(IB2) * TMPX
                M(1,2) = M(1,2) - KX * Weight(IB2) * TMPZ
                M(2,1) = M(2,1) - KZ * Weight(IB2) * TMPX
                M(2,2) = M(2,2) - KZ * Weight(IB2) * TMPZ
            Endif
          Enddo
        Enddo
        Enddo
        Enddo
        !Write(*,*) M(1,1),M(1,2), M(2,1), M(2,2)

      Call Matrix2Dinv(M)
      
      End Subroutine
